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The theoretical description of transport in a wide class of novel materials is based upon quan- 
tum percolation and related random resistor network (RRN) models. We examine the localization 
properties of electronic states of diverse two-dimensional quantum percolation models using exact 
diagonalization in combination with kernel polynomial expansion techniques. Employing the lo- 
cal distribution approach we determine the arithmetically and geometrically averaged densities of 
states in order to distinguish extended, current carrying states from localized ones. To get further 
insight into the nature of eigenstates of RRN models we analyze the probability distribution of the 
local density of states in the whole parameter and energy range. For a recently proposed RRN 
representation of graphene sheets we discuss leakage effects. 
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I. INTRODUCTION 

Disorder is an intrinsic feature of many solid state ma- 
terials. The spatial inhomogeneity of a sample strongly 
affects the transport properties, particularly in low- 
dimensional systems where Anderson localization effects 
play an important role 1 . Attempts to describe quan- 
tum transport in disordered media usually rely on free 
electron models with random links between lattice sites 
and/or varying on-site potentials. Then transport is re- 
lated to percolating current patterns in a kind of random 
resistor network (RRN). RRN models apply to both clas- 
sical and quantum percolation problems. For classical 
percolation the existence of a sample-spanning cluster 
ensures finite conductivity above a certain percolation 
threshold of accessible sites p c . In the quantum case, the 
appearance of current carrying states strongly depends 
on the (relative) importance of tunneling, scattering and 
interference processes. For instance, strong scattering at 
the irregular boundaries of a conducting cluster may lead 
to a localization of the charge carrier even for p > p c . 

Over the past years percolation approaches have been 
employed in many circumstances, e.g., in order to model 
the transport properties of doped semiconductors 2 - and 
granular metals 3 , the dynamics of atomic Fermi-Bose 
mixtures 4 , the wave propagation through binary inho- 
mogeneous media 5 , as well as the formation of novel 
two-dimensional spin liquids 6 . Another focal point is 
the metal-insulator transition, e.g., in 2D n-GaAs het- 
erostructuresi and colossal magnetoresistive manganite 
films 8 , or rather the superconductor-insulator and (in- 
teger) quantum Hall transitions^^. Quite recently, the 
problem of disorder in systems with Dirac fermions has 
been studied in the context of dirty superconductors^ 
and two-dimensional (2D) graphen o 12 i 13 i 14 i 15 . 

In mono- and bilayer graphene, and graphene-based 
field effect transistors, strong fluctuations of the lo- 
cal charge density and, thus, the local conductiv- 
ity have been reported near the so-called "neutrality 



point 



» 16.17.18.19 



where the conductivity reaches its min- 



imum. The mesoscopic regions of different charge car- 



rier density may be caused by inhomogeneities in the 
substrate, or non-perfect stacking 20 i 21 i 22 i 23 . In order to 
model the minimal conductivity in graphene a RRN rep- 
resentation of a graphene sheet has been proposed by 
Cheianov et al. 2A , where random links between electron 
and hole "puddles" (corresponding to lattice sites) are 
assumed to determine the observable conductivity rather 
than the local conductivity of a puddle. Such a RRN 
formulated on a square lattice is closely related to a 2D 
quantum percolation percolation model with additional 
finite "leakage" between all lattice sites. 

Motivated by this situation, in this paper we per- 
form an in-depth numerical study of generalized 2D per- 
colative RRN models. In particular we analyze how 
the leakage and connectivity rules influence the nature 
of the single-particle states. Since in two dimensions 
the problem of quantum percolation is still discussed 
controversially — especially with respect to the existence 
of a quan t um per colation threshold p c < p q < 1, see, e.g., 
Refs^^^^^^^^l^^M—we rely on un- 
biased numerical techniques, which take quantum effects 
fully into account. To this end we employ the so-called 
local distribution (LD) approac h 37 ! 38 . This technique, 
based on the determination of the distribution of the local 
density of states (LDOS) for all energies, has been previ- 
ously applied to tackle localization phenomena in various 
disordered systems with great succes o 39 i 40 i 41 i 42 i 43 i 44 i 45 i 46 . 

The paper is organized as follows. Sect. II introduces 
the RRN model and briefly outlines the LD approach. 
In Sect. Ill, we examine the localization properties of 
the eigenstates of RRN models by calculating the LDOS 
distribution, as well as the average and typical density of 
states. Sect. IV contains our conclusions. 



II. MODEL AND METHODS 

Let us consider a 2D square lattice with TV sites on 
two sub-lattices a and (3, which, e.g., might represent 
regions of different charge carrier concentrations. The 
sub-lattices are linked by connection rules characteris- 



tic of specific materials. Generation rule A corresponds 
to a checkerboard- like structure (see upper left panel of 
Fig. [I]). Regions (sites) of the sub-lattices are randomly 
connected with each other (lower panels of Fig. [1]) ac- 
cording to the connection rules displayed in the upper 
right panel of Fig. [H The hopping probability between 
suchlike connected sites is assumed to be much higher 
than for hopping events to nearest neighbors. The latter 
ones are reduced by the leakage A < 1. 

For the more abstract case of generation rule B (with 
the above connection rule), additional random bonds 
connect a to (3 sites, favoring the formation of quasi one- 
dimensional zig-zag chains. We still have a weak leak- 
age between all neighboring sites. This generation rule 
can be interpreted as an attempt to incorporate a spatial 
anisotropy into the model. 

In both cases, the corresponding RRN can be described 
by the Hamiltonian 



H 



+ Yj (W+T C *+- + ( X - Vi)4 c i+S 

ie/3 
+ X Y1 (4 c *+-4 c m) 



(1) 



H.c. . 



Here c] (q) creates (annihilates) an electron at site i and 
the arrows denote the nearest neighbors site of i in the 
corresponding direction. The discrete random variables 
r]i G {0, 1} determine which diagonal in the plaquette is 
connected (cf. upper right panel in Fig.[T]). By fixing the 
expectation value of the {?7i}-distribution, p — (rji), we 
can control the sizes of the regions of connected lattice 
sites. For p = 0.5 the RRN model (pQ) does not depend on 
the generation rules A or B because of the symmetry of 
the bimodal {?7i}-distribution. In the limit of vanishing 
leakage we obtain the standard 2D quantum percolation 
model using generation rule A. 

To characterize the transport behavior of our model we 
adopt the local distribution approach (for a detailed de- 
scription of this technique, we refer the reader to Ref. |47| 
and references therein). The main idea is the following: 
We calculate the LDOS for all lattice sites i, 



N 



Pi (E) = J2\(i\n)\ 2 S(E-E n ). 



(2) 



This quantity depends on the energy and varies from site 
to site. Moreover the pi are specific for a particular sam- 
ple (realization {rji}). Calculating and recording pi for 
many sites and realizations we obtain the probability dis- 
tribution f[pi] and distribution function 



F[Pi 



MWi 



(3) 



of {pi}. Both f[pi] and F[pi] are self- averaging. That is, 
in the thermodynamic limit, they are independent of the 
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FIG. 1: (Color online) Upper left panel: Generation rules A 
and B for sub-lattices a and f3. Upper right panel: Connec- 
tion rule of neighboring sites on different sub-lattices. Lower 
panels: Three particular cluster realizations representing the 
above rules on a 20 x 20 lattice. From left to right: p = 0.50, 
p = 0.90 (rule A), p = 0.90 (rule B). Note that rules A and B 
become equivalent for p = 0.50. At the red dashed lines the 
system is interconnected by periodic boundary conditions. 



actual realization {rji} and the chosen sites {z}, but will 
be solely determined by the (global) model parameters 
p, A, and the generation rule. This restores translational 
invar iance on the level of distributions. 

A well-established criterion for localization then can 
be deduced from the shape of the distribution function. 
Since for an extended (current carrying) state the am- 
plitude of the wave function is more or less uniform, 
F[pi] steeply rises in the vicinity of the mean value of 
the LDOS, p me = (pi). For localized states, on the other 
hand, the LDOS strongly fluctuates throughout the lat- 
tice leading to a rather gradual increase of F\pi\ as a 
function of the magnitude of pi . To capture this different 
behavior quantitatively, we may compare the arithmetic 
mean of the LDOS, p me , to its the geometric mean, the 



so-called "typical" DOS, p ty = e 



= P (Hpi)) 



, see, e.g., Ref. I3S 



For extended states p me and p ty are of the same order, 
whereas p ty vanishes for a localized state (or, in a finite 
system, is at least considerably suppressed). Of course, a 
reliable distinction between extended and localized states 
requires the consideration of different system sizes. For 
extended states F[pi] is independent of N. By contrast, 
for localized states F[pi] shifts toward lower values as N 
increases (causing the reduced value of p ty ). 

The LDOS can be calculated very efficiently by means 
of the kernel polynomial method (KPM) 48 . Within this 
technique, the spectrum of H is expanded into a finite se- 
ries of Chebyshev polynomials with additional damping 
factors. This approximation can be viewed as a convolu- 
tion of the spectrum with the Jackson kernel, an almost 
Gaussian of width a. As a depends on the order of the 
Chebyshev series as well as on the energy, we have to 
adapt the expansion order to ensure an uniform resolu- 



tion (constant a) for the whole spectrum^. 



III. RESULTS 

In Fig. [2] we compare the mean and typical DOS in de- 
pendence on A and p for generation rules A and B. While 
the 2D percolation model (A = 0, rule A) exhibits sym- 
metric DOS spectra, the inclusion of next-nearest neigh- 
bor hopping causes an asymmetry that grows with in- 
creasing A. Moreover, if p > 0.5, the results obtained 
for the two generation rules differ significantly. For the 
A-case the mean DOS resembles the 2D DOS. Genera- 
tion rule B leads to a spectrum evocative of the DOS 
in ID. In both cases a multitude of spikes exists which 
we can attribute to localized states on "isolated" islands. 
This feature is well known from the binary alloy model 49 . 
Increasing the leakage broadens the peaks and reduces 
their abundance, because the "isolation" of the islands is 
weakened. Furthermore, the leakage shifts the "special 
energies" at which these peaks appear toward lower val- 
ues as compared to those in the standard 2D percolation 
model^ 9 -. 

The typical DOS underlines the prominent role of the 
leakage. In Fig. [2] (a) a vanishing p ty suggests that all 
states are localized. As long as the leakage is very small, 
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FIG. 2: (Color online) Mean (upper solid lines) and typical 
(lower dashed lines) DOS. Results are shown for N = 400 2 
and Na = 45 for different p, A and generation rules A and B. 
To illustrate the finite-size dependence of p ty the results for a 
N = 800 2 lattice are included (long dashed lines). The data 
are based on 10 5 (10 4 ) individual LDOS spectra for N = 400 2 
(800 2 ). Dotted vertical lines indicate the energies for which 
the distribution functions of the LDOS are given in Fig. [3] 



the typical DOS remains small as well throughout the 
band. This also holds for larger p and both generation 
rules. Sizeable values of p ty appear at larger values of A 
only (compare the data for A = 0.01 and 0.1). As for the 
mean DOS, we observe a pronounced asymmetry of pt y in 
the electronic band. This favors a finite typical DOS, i.e. 
extended states, for E > 0. However, this has to be taken 
with caution. If we increase the system size together with 
the resolution of the KPM in such a way, that the number 
of states within the width of the Jackson kernel, Na, 
is kept constant, the picture changes. While for truly 
extended states, p ty should be independent of the system 
size (cf. Fig. [2] (d)), in Figs. [2] (b) and (f), the typical 
DOS decreases with increasing system size. This points 
towards a large localization length. While for small and 
moderate N the state still spans the entire lattice, the 
system size exceeds the localization length for larger N. 
Therefore, on a considerable number of sites the LDOS 
is very small, leading to a reduced value of p ty . Then 
the question arises, why in panel Fig. [2] (d) the typical 
DOS is reduced as compared to p me , even though this 
ratio is independent of N. This can be easily understood 
if we consider a perfectly ordered system, i.e. p = 1. 
There the sites which do not belong to the spanning sub- 
lattice are also taken into account, and amount to half 
of the lattice sites. In absence of leakage, at those sites 
the LDOS vanishes and probing such a site completely 
suppresses p ty . For finite A, however, these sites have 
small amplitudes and in this manner p ty is reduced as 
compared to p me but finite. These arguments also hold 
for panel Fig.[2](b), but there the low leakage suppresses 
the typical DOS already almost completely, except for 
E = 0. 

As mentioned above, the comparison of p ty and p me 
may only serve as a first indication of localization. A 
careful finite-size scaling of the full probability function 
(or the distribution function) of the LDOS allows for a 
more reliable distinction between extended and localized 
states. Thereto, in Fig. [3] we show for the same sets 
of parameters the distribution function of the LDOS at 
various characteristic energies. In the panels of Fig. [2] 
these energies are indicated by vertical dotted lines. 

Here one crucial point concerns the normalization of 
pi to p me . For a completely extended state, with uni- 
form amplitudes, the distribution function would be a 
step function. Due to the used scaling, this step would 
be located at Pi/p m e = 1, irrespective of the system size. 
The random structure of the underlying lattice will dis- 
tort this perfect jump-like behavior of F[pi/p me \. But 
nevertheless extended states are characterized by an N- 
independent distribution function [cf. Fig. [3] (d)]. In 
contrast, for localized states, F[pi/p me ] constantly shifts 
toward lower values as TV increases [cf. Figs. [3] (a), (e)]. 
Depending on the localization length, this effect is more 
or less pronounced. Thus, a wide range of system sizes 
is necessary in order to discriminate localized from ex- 
tended states by means of finite-size scaling [cf. Figs. [3] 
(b) and (f)]. 




FIG. 3: (Color online) Distribution function of the normalized 
LDOS for the energies indicated by vertical lines in Fig. [2] 
The statistics is based on 10 7 , 10 6 , 10 5 , 10 4 LDOS values for 

AT: 



100 2 , 200 2 , 400 2 , 800 2 and resolution Na ■■ 



45. 



A particular interesting shape of F[pi/p me ] appears in 
Fig. 03(c), for a state in the band center. From the multi- 
step structure, we may deduce that the probability distri- 
bution is mainly concentrated around three values. The 
largest of those values is independent on A/", while the 
others show the above discussed variation. This behav- 
ior can be explained by considering again the completely 
ordered case. In absence of leakage the E = eigenstate 
is highly degenerate: we have N/2 completely localized 
states, one on each isolated site. The other half of the 
eigenstates are extended in the perfect lattice, providing 
energies in the whole band. As the LDOS probes the 
complete eigenspace and not just the amplitude of one 
particular eigenstate, at E — we obtain the same value 
of pi/pme for each isolated site. Introducing a small leak- 
age kills the high degeneracy in principle. Due to the 
finite energy resolution of the KPM, however, still many 
of these states contribute to the LDOS at E = 0. In 
the presence of imperfections (p < 1), some of the for- 
mer isolated sites will now be connected to form larger 
islands. Then less than N/2 states will contribute to the 
"isolated islands" peak at high Pi/p me [but still around 
38% in Fig. [3] (c)]. The second step, having a weight of 
about 50%, is due to sites on the highly connected major- 
ity sub-lattice, on which the LDOS is reduced but finite 
due to the leakage. The remaining lowest step originates 
from more complicated islands of several sites. On those 
islands (2, 3, . . . sites), the E = eigenstates have vanish- 
ing amplitudes on some sites. Due to the leakage these 
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FIG. 4: (Color online) Normalized occupation probability 
N\{i\n)\ 2 of several characteristic eigenstates \n) on a N = 
128 2 lattice. The generation rules, doping and leakage as well 
as the chosen energies E n correspond to the ones for which 
the distribution functions of the LDOS has been presented 
in Fig. [3] (same order of panels). In the center column the 
distribution function of the normalized occupation probabil- 
ity is shown for these states. Furthermore, for each state an 
enlargement of a characteristic region together with the local 
structure of the present links is shown. 



sites again acquire a finite value of pi . 

To get additional insight into the spatial structure of 
the eigenstates, we investigate a smaller 128 x 128 sys- 
tem using exact diagonalization techniques. Figured! vi- 
sualizes characteristic eigenstates for the parameters dis- 
cussed in Figs. [2]and [3l Most notably is the pronounced 
leakage-dependence of the results. While we see clear lo- 
calization for A = 0.01 [Fig. [4] (a), (e)], in the high- leakage 
case all states span the entire lattice. The ID structure 
for rule B and large p, which we already found in the 
DOS in Fig. [2] (e), is also predominant in Fig. [4] (e). In 
Fig. 21(c), notably in the magnifying inset, we see our as- 
sumption confirmed, that sites with large amplitudes are 
almost exclusively isolated, with the exception of some 
larger islands. Having the sharp step in mind, which oc- 
curs in Fig. [3] (c) at large Pi/p me , the alert reader may 
wonder about the different amplitudes in Fig. [4] (c) and 
the rather smooth increase of F[AT|(i|n)| 2 ] in the mid- 
dle column. This discrepancy is again due to the high 
degeneracy of the E = state and the fact, that the 
LDOS takes into account the whole E = eigenspace. 
In contrast, for A/"|(z|n)| 2 , we randomly pick one particu- 
lar eigenstate out of this subspace. For the other E = 
eigenstates, F[N\(i\n)\ 2 ] looks similar, only the location 
of the sites with maximum amplitudes differ. Clearly, the 



physically relevant quantity is the LDOS, as the choice of 
eigenvectors which span the eigenspace is arbitrary up to 
any linear combinations of them. Summing up all eigen- 
states to E = for a fixed lattice site we indeed get the 
same amplitude on each isolated site, which brings us 
back to the findings of Fig. [3] (c). Finally we consider 
the case of higher leakage (A = 0.1), where all states 
look rather similar. Here the amplitudes fluctuate over a 
smaller range, and the fluctuations occur on very short 
spatial scales. This is most pronounced in Figs. [H (d), 
(f) where no global structures can be distinguished, in 
accordance with the notion of extended states. By con- 
trast, although the state in Fig. [H (b) spans the entire 
lattice, we observe sizeable regions with higher and lower 
amplitudes than the mean value. Presumably, these in- 
homogeneities are even more pronounced for larger sys- 
tems which, however, (as yet) are not accessible by our 
exact diagonalization studies. In any case, the shifting 
of the LDOS distribution function in Fig. [3] (b) suggests 
this state to be localized on large length scales. 



bining exact diagonalization, Chebyshev expansion and 
local distribution approaches in calculating the local den- 
sity of states, we determine — after a careful finite-size 
scaling — the localization properties of the single-particle 
eigenstates. We found that current carrying states ex- 
ist, and that the appearrance of diffusion is mainly trig- 
gered by the amount of leakage. In contrast to previous 
work we analyzed the nature of the single-particle states 
in the whole energy band and not just the vicinity of 
the neutrality point in the band center E = 0. A ten- 
dency toward extended states is observed for E > 0. In 
view of the simplicity of the considered RRN model, our 
data are certainly not yet suited to be compared against 
real transport data for, e.g., undoped graphene mono- 
layers. Nevertheless, using unbiased numerics on a mi- 
croscopic level, we fully account for quantum effects in 
a RRN model originally designed to describe transport 
in a graphene-based field effect transistor, and received 
results that support a minimal conductivity in graphene. 



IV. SUMMARY 
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